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ABSTRACT 


The  structure  and  behavior  of  barotropi cal ly  unstable 
and  stable  wave  disturbances  in  the  vicinity  of  a zonally 
varying  easterly  jet  are  studied  numerically  with  a 
linearized  barotropic  vorticity  equation  on  a 8-plane.  The 
easterly  jet  is  approximated  by  a Bickley  jet  with  a slow 
zonal  variation.  The  numerical  results  are  also  compared 
with  a simple  mechanistic  analytical  model  using  the  local 
phase  speed  and  growth  rate  concepts.  The  results  are 
grossly  similar  in  several  respects  to  that  expected  from 
the  parallel  flow  theory  of  barotropic  instability,  however, 
the  resultant  structure  of  the  waves  causes  a spatial  growth 
rate  greater  than  predicted  by  the  local  growth  rates 
computed  wtth  a parallel  flow  model.  In  the  stable  region, 
the  structure  leads  to  strong  dynamic  damping.  When  a 
uniform  advective  velocity  is  added  to  a variable  mean  flow, 
the  differences  between  the  behavior  of  the  computed  waves 
and  that  implied  by  the  parallel  flow  theory  are  somewhat 
reduced.  The  waves  remove  kinetic  energy  from  the  mean 
flow  and  most  of  this  energy  is  removed  on  the  downwind  side 
of  the  jet.  The  computed  structure  and  behavior  of  the  waves 
have  a number  of  features  that  resemble  those  observed  in 
the  vicinity  of  the  upper  troposphere  easterly  jet  during 
the  summer  monsoon  . 
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I.  INTRODUCTION 


In  the  Northern  Hemisphere  summer,  a moderately  strong 
easterly  jet  exists  south  of  the  Tibetan  high  near  the  200  mb 
level  (Kri shnamurti  , 1971a,  1971b).  Synoptic-scale  moving 
dfsturbances  occur  at  the  level  of  the  jet,  and  it  appears 
that  these  disturbances  arise  from  barotropic  instability  of 
the  mean  flow.  The  jet  contains  regions  of  large  vorticity 
gradients  where  the  necessary  condition  for  barotropic 
instability  is  sometimes  locally  satisfied.  If  the  observed 
disturbances  were  the  result  of  barotropic  instability,  they 
would  extract  energy  from  the  mean  zonal  flow  and  the  plane- 
tary-scale waves,  since  the  latter  combine  with  the  zonal 
flow  to  give  the  large  vorticity  gradients  south  of  the 
Tibetan  high.  In  fact,  Kanamitsu  e_t  aj_.  ( 1 972  ) have  shown 
that  wavenumbers  6-8  in  the  wind  spectrum  in  the  region 
between  1 5 S and  15N  receive  energy  through  barotropic 
interaction  with  the  zonal  and  wavenumber  1 flow. 

Furthermore,  Krisrinamurti  (1971a)  studied  the  wavenumber 
spectra  of  the  meridional  wind  component  for  selected 
tropical  latitudes  and  he  found  a peak  in  the  spectra  near 
wavenumbers  6-8  at  latitudes  near  the  easterly  jet. 

Colton  (1973)  studied  barotropic  interactions  between 
quas i -s tati onary  long  waves  and  transient  synoptic  waves 
using  a semi -spectral  numerical  model.  His  long  waves  were 
forced  with  a specified  divergence  field  following  the 


r 


diagnostic  model  of  Holton  and  Colton  (1972).  Results  from 
a long  term  integration  of  his  model  reproduced  several 
features  of  the  upper  tropospheric  general  circulation.  In 
particular,  his  model  simulated  the  disturbances  observed  by 
Krishnamurti  (1971a,  1971b)  in  the  vicinity  of  the  easterly 
jet.  His  model  showed  that  wave  disturbances  entering 
upstream  of  the  easterly  jet  regime  grow  and  increase  their 
speed  and  wavelength  as  they  approach  the  longitude  of  the 
jet  velocity  maximum,  achieve  their  maximum  amplitude  a 
considerable  distance  downstream  of  the  jet  maximum,  and 
eventually  decay  as  the  disturbances  move  out  of  the  easterly 
jet  regime.  Fur the  rmore , the  high  and  low  pressure  centers 
of  these  transient  wave  disturbances  occur  on  the  wings  of 
the  central  axis  of  the  jet  stream.  Colton  concluded  that 
the  transient  disturbances  are  due  to  scale  interactions 
Vnvolvtng  short-term  barotropic  instability. 

In  both  Kri s hnamurti  1 s observational  and  Colton's 
numerical  studies,  the  zonal  variation  of  the  jet  apparently 
has  significant  effects  on  the  dynamic  behavior  of  the 
transient  disturbances.  This  study  is  an  attempt  to  under- 
stand better  such  effects.  This  is  a linear  study,  therefore, 
the  results  can  only  be  suggestive  of  the  complicated  nonlinear 
behavior  of  the  real  atmosphere.  However,  since  the  ampli- 
tudes of  synoptic  waves  are  generally  of  smaller  magnitude 
in  the  tropics  than  in  middle  and  higher  latitudes,  a linear 
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stability  approach  should  be  a reasonanle  approx i mati on  to 
the  real  atmosDhere. 
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This  study  is  an  extension  of  the  linear  barotrooic 
instability  theory  developed  over  the  years  by  various 
investigators.  The  extension  of  this  study  is  the  zonal 
variation  of  the  basic  flow.  Rayleigh  (1890,  1913)  developed 
the  concept  of  hydrodynamic  instability  and  Kuo  (1949) 
extended  this  concept  to  a rotating  atmosphere  by  including 
the  beta  term.  Since  then,  many  investigators  have  studied 
the  stability  of  barotropic  zonal  flows.  Kuo  (1949,  1951), 
Lipps  (1962),  Yanai  and  Nitta  (1968)  investigated  barotropic 
instability  of  parallel  symmetric  westerly  jet  flows.  Nitta 
and  Yanai  (1969),  Yamasaki  and  Wada  (1972),  and  Kuo  (1973) 
extended  this  theory  to  parallel  symmetric  easterly  jet  flows. 
Lipps  (1965,  1970)  and  Kuo  (1973)  further  extended  the  concept 
to  parallel  asymmetric  barotropic  zonal  flows.  In  general, 
when  the  necessary  criterion  for  linear  barotropic  instability 
is  satisfied,  normal  mode  solutions  for  barotropic  unstable 
waves  have  certain  properties:  1)  they  exist  within  a range  of 
intermediate  wavelengths,  2)  they  have  a latitudinal  tilt 
opposite  to  the  shear  of  the  barotropic  zonal  wind  and  this 
tilt  creates  momentum  fluxes  which  redistribute  the  kinetic 
energy  from  the  mean  zonal  flow  to  the  disturbance,  3)  their 
phase  speeds  in  a westerly  jet  are  less  then  the  speed  of  the 
maximum  westerly  wind,  and  4)  their  phase  speed  in  an 


easterly  jet  may  be  greater  than  the  maximum  easterly  wind 
(Pedlosky,  1964;  Yamasaki  and  Wada,  1972).  In  this  case, 
however,  the  growth  rates  are  generally  very  small. 

Nitta  and  Yanai  (1969)  modified  the  concept  that 
unstable  solutions  exist  only  within  an  intermediate  range 
of  wavelengths.  For  an  easterly  sinusoidal  jet  flow,  they 
found  a distinct  short  wave  cut-off  but  no  apparent  long  wave 
cut-off  for  instability.  Yamasaki  and  Wada  (1972)  modified 
this  by  showing  that  the  long  wave  cut-off  is  dependent  on 
the  strength  and  sharpness  of  this  sinusoidal  easterly  jet. 
For  very  strong  and  sharp  velocity  profiles,  the  long  wave 
cut-off  approaches  infinity.  For  relatively  weak  but  still 
unstable  profiles,  a finite  long  wave  cut-off  exists  for 
this  jet.  Kuo  (1973)  referred  to  this  long  wave  cut-off 
region  as  the  modified  Rossby  regime  in  his  numerical  study 
for  an  easterly  Bickley  jet. 

All  these  stability  studies  have  used  parallel  baro- 
tropic  flows,  i . e . , there  is  no  zonal  variation  in  the  basic 
flow.  Lorenz  (1972),  however,  investigated  primarily  by 
analytical  means  the  barotropic  instability  of  a flow  pattern 
which  varies  with  longitude.  The  basic  flow  is  a neutral 
Rossby  wave  superposed  on  a uniform  westerly  flow.  Essen- 
tially, this  flow  depicts  the  progression  of  large  scale 
waves  embedded  in  a westerly  current.  Zonal  flows  of  mid- 
latitudes are  generally  considered  to  be  barotropi ca 1 ly 
stable,  but  Lorenz  showed  that  a uniform  zonal  flow  together 
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with  superposed  neutral  Rossby  waves  may  be  unstable  with 
respect  to  further  perturbations.  Lorenz  found  that  the 
growth  rate  of  the  perturbations  is  comparable  to  the  growth 
rate  of  the  errors  of  large  numerical  atmospheric  models. 

Based  on  this,  Lorenz  suggested  that  barotropic  instability 
may  be  partly  responsible  for  the  unpredictability  of  the 
real  atmosphere. 

In  addition  to  the  regular  or  normal  mode  ei gensol uti ons  , 
there  are  "singular"  or  "continuum"  mode  solutions  which 
have  been  discussed  by  Rayleigh  (1913).,  Case  ( 1 960),  Pedlosky 
(1964),  and  Yanai  and  Nitta  (1968).  These  "singular" 
solutions  correspond  to  continuous  eigenvalues  of  phase 
velocity  (c)  which  are  equal  to  the  basic  flow,  IT(y),  some- 
where in  the  zonal  current.  Since  these  "singular"  solutions 
are  continuous,  there  is  an  infinite  number  of  solutions. 

Case  showed  that  these  continuum  modes  are  needed  to  form  a 
complete  set  of  solutions.  He  also  showed  that  the  disturb- 
ance formed  by  the  sum  of  the  continuum  modes  has  a 
y-structure  tilt  in  the  same  sense  as  the  basic  flow  shear 
and  that  it  usually  decays  as  1/t  or  faster,  where  t repre- 
sents time. 

It  is  clear  that  in  the  real  atmosphere,  "mean  flows" 
vary  both  in  space  and  time.  They  are  neither  purely 
barotropic  nor  purely  baroclinic,  and  in  including  both 
effects  in  a linear  stability  study  is  a very  difficult  task. 
Most  of  the  investigators  have  examined  the  dynamic 

II 
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instability  problem  by  studying  each  effect  separately. 

This  study,  however,  focuses  only  on  baro tropic  instability. 
We  hope  that  a better  understanding  of  the  kind  of  upper 
tropospheric  waves  studied  by  Krishnamurti  (1971a,  1971b) 
and  Colton  (1973)  near  the  easterly  jet  can  be  achieved  by 
studying  the  behavior  of  barotropic  waves  in  a region  of  a 
variable  mean  easterly  wind.  We  are  especially  interested 
in  the  situation  where  the  waves  move  into  and  out  of  baro- 
tropically  unstable  regions.  Therefore,  the  objective  of 
this  work  is  to  determine  the  dynamic  stability  behavior, 
the  structure  and  the  energetics  of  barotropic  waves  propa- 
gating in  a zonally  varying  mean  wind. 

With  a very  simple  analytical  model,  Pedlosky  (1976) 
studied  the  finite  amplitude  dynamics  in  a zonally  varying 
baroclinic  current.  The  mean  flow  of  his  model  changes 
abruptly  from  a weakly  unstable  regime  with  a constant  growth 
rate  to  a stable  regime  downstream.  He  found  that  disturb- 
ances may  propagate  into  stable  regions  substantially 
undiminished,  retaining  a considerable  memory  of  its  history 
in  a locally  unstable  region. 

The  mean  flow  in  our  barotropic  study,  however,  has  a 
full  downstream  variation  such  that  the  stability  of  the 
mean  flow  ranges  from  locally  strong  instability  to  locally 
strong  stability. 


1 7 
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Our  approach  is  as  follows:  firstly,  we  develop  a 
numerical  model  based  on  the  linear  non-di vergent  barotropic 
vorticity  equation  on  a beta-plane  incorporating  the  dynamics 
of  transient  barotropic  waves  in  a region  of  variable  mean 
wind.  The  men  *** i ^ d is  an  easterly  hyperbolic  secant- 

squarred  (Bick  .d  the  mean  meridional  wind  is  derived 

in  such  a way  that  che  mean  flow  is  non-di vergent . Secondly, 
using  th’s  numerical  model,  we  conduct  selected  case  studies 
to  determine  fundamental  dynamic  stability  and  energetic 
properties  of  these  barotropic  waves.  Thirdly,  we  compare 
the  results  obtained  from  the  numerical  model  with  the 
parallel  flow  theory.  This  is  done  by  constructing  a simple 
me chan i Stic  analytical  model  which  incorporates  the  local 
stability  concept  of  the  parallel  flow  theory.  This  compari- 
son gives  us  further  insights  on  the  effect  of  the  variable 
mean  wind. 

The  domain  of  the  numerical  model  consists  of  an  open 
channel  with  rigid  walls  at  the  north  and  south  boundaries. 
The  basic  flow  is  a slowly  varying  easterly  Bickley  jet. 

A periodic  forcing  is  applied  on  the  inflow  or  eastern 
boundary  such  that  periodic  perturbations  are  generated  from 
this  boundary  into  the  channel  flow  which  represents  waves 
moving  into  the  region  from  the  east.  The  western  boundary 
condition  allows  the  waves  to  move  out  of  the  channel.  As 
the  waves  move  through  the  region,  they  grow  or  decay  in 
relation  to  the  local  stability  properties  of  the  mean  flow, 
whereas  at  each  point  the  fields  vary  periodically. 
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II  . NUMERICAL  MODEL 
A.  GENERAL  FORMULATION 

The  governing  equations  are  the  barotropic  vorticity 
equation  using  a beta  (a)  plane  approximation: 

fl  * ulf  + vly  + SV  - Q - DfS  (2.1) 

and  the  non-di vergent  continuity  equation 


3 U 3 V 

3 x 3y 


0 , 


(2.2) 


where 


_ 3 V 3 U 
5 = 3X  * 3y 


(2.3) 


Here  3 is  the  north-south  gradient  of  the  earth's  vorticity 
and  is  given  by  its  value  at  10°  latitude.  Q is  a forcing 
function  representing  non-barotropi  c effects  which  is 
required  to  maintain  the  vorticity  field,  is  a frictional 
coefficient.  Since  the  flow  is  two-di mens i onal  and  non- 
divergent,  a streamfunction  (^ ) is  defined  such  that 


IA  - 11 

3y  ’ 3 x * 


(2.4) 


thus 


t . (2.5) 

Here  u and  v are  the  velocities  of  the  flow  in  the  x and  y 
direction,  respectively.  Equation  (2.1)  now  becomes 


J_  72j,  11  + 11  + aH 

3 1 ^ ~ 3y  3 X 3 X 3y  a 3 X 


Q 


Df7% 


(2.6) 
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The  governing  equations  are  linearized  by  separating 
quantities  into  the  following  form: 

<*(x,y,t)  = $(x,y)  + p ' ( x , y , t ) . (2.7) 


Here  the  bars  are  the  basic  state  variables  and  the  primes 
are  the  perturbation  quantities.  The  mean  s treamf uncti on 
i£(x,y)  satisfies  the  mean  vorticity  equation  when  a non- 
barotropic  source  term  is  added.  The  basic  state  equation 
is  therefore 


. li  3 7 2p  + 3$.  3 + 3j£.  _ ; 

ay  3 X 3 X 3y  3x  w 


(2.8) 


which  is  non-linear.  Q is  needed  to  maintain  the  basic 
vorticity  field  in  steady  state.  The  resultant  linearized 
perturbation  equation  is: 
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3 t 
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3 X 
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= - D^7  p 

3?2" 
3 x 


(2.9) 


The  advective  terms  in  (2.9)  can  be  written  in  Jacobian  form 
such  that 


= - J(*,vV)-  JU'.v2’})  - sf-^—  - Df72p'  , (2.10) 

where  (2.10)  is  a Poisson  equation  with  the  tendency  of  the 
perturbation  stream-function  as  the  dependent  variable.  The 
right  side  of  (2.10)  is  the  forcing  or  non-homogeneous  oart 
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of  the  Poisson  problem  for  3^'/3t.  The  north-south  boundaries 
are  ri gi  d walls  where 


- 0 , 


(2.11) 


and 


(2.12) 


The  inflow  or  eastern  boundary  conditions  are  specified  such 
that 


♦'(O.y.t)  = A (y ) Sin(ut)  + 3(y)  cos(ut),  (2.13) 

and 


(O.y.t)  = A(y)ucos(ut)  - s (y )us 1 n (ut ) , (2.14) 

where  a is  the  specified  frequency  of  the  forcing.  A(y) 
and  B(.y)  are  the  y-structure  coefficients  of  the  forcing 
which  are  described  in  Chapter  IV.  The  purpose  of  this 
forcing  is  to  introduce  waves  into  the  region  from  the  east. 
The  other  condition  on  the  inflow  boundary  is 

2 

c'(0,y,t)  = vV  « - kV  + . (2.15) 

3y 

This  expression  follows  because  it  is  expected  that  $ 1 will 

i k x 

have  a spatial  variation  of  the  form  e . The  determination 
of  k is  discussed  later.  This  boundary  forcing  is  expected 
to  lead  to  periodic  wave  disturbances  throughout  the  domain 
after  a certain  time  period  of  numerical  integration.  It 
turns  out  that  the  outflow  boundary  conditions  are  crucial 
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for  attaining  this  periodic  state,  because  the  wave  disturb- 
ances must  be  able  to  propagate  freely  out  of  this  boundary. 
Therefore,  a Sommerfel d radiation  condition  (Pearson,  1974) 
for  both  the  tendency  and  the  vorticity  is  used  to  approxi- 
mate this  mechanism  on  the  outflow  boundary: 

ft  (If’  ■ - cr.  tx  <ff).  ■"«* 


ill  . at' 

3t  rQ  3 X 


(2.17) 


Here  c is  a specified  constant  phase  velocity.  Pearson 
' o 

(1974)  showed  that  the  Sommerfeld  radiation  condition  is  a 
consistent  boundary  condition  for  numerical  models  of  initial 
value  systems  admitting  dispersive  waves.  The  fundamental 
problem  with  most  outflow  boundary  conditions  is  the  reflec- 
tion of  incident  waves  from  a boundary  back  into  the  interior 
region.  This  is  usually  disastrous  for  numerical  models 

admitting  dispersive  waves.  Pearson  showed  that,  if  c^  is 

o 

chosen  judiciously,  the  longwaves  move  smoothly  through  the 
boundary.  For  short  wavelengths,  however,  an  area  with  a 
large  coefficient  of  viscosity  near  the  boundary  is  often 
required  to  help  control  the  reflection  problem.  According 
to  Pearson,  the  amount  of  damping  is  proportional  to  the 
wavenumber.  A Matsuno  finite  difference  scheme  is  used  in 
our  numerical  model.  This  scheme  has  a tendency  to  damp  the 
short  waves  thus  we  need  only  to  use  a relatively  small 
coefficient  of  friction  to  control  wave  reflection  at  the 


outflow  boundary.  A friction  coefficient  of  0^  = .15x10“ 
sec"^  is  used  which  is  equivalent  to  an  e-folding  decay  time 
of  approximately  7.7  days. 


B.  ZONALLY  VARYING  BASIC  FLOW 

The  basic  velocity  field  is  an  easterly  Bickley  jet 
defined  by 

u(x,y)  = - U ( x ) sech2(j^y)  - UQ  = - |y  . (2.18) 

Here  d(x)  is  a characteristic  length  scale  of  the  jet  and  is 
related  to  the  half  width  d(x)  (Kuo,  1973)  by 

d ( x ) = - 1.76  d (x ) . (2.19) 


U(x)  is  the  velocity  of  the  jet  at  y=0  and  UQ  is  a constant 
velocity.  The  basic  flow  streamfunction  ( H ) is  specified 
to  be  constant  at  both  the  southern  (y=-D)  and  northern  (y=D) 
boundaries  as  follows: 


’J(-D)  = C,) 

, where  (2.20) 

5(D)  = c2  j 

(2.21) 

y y 

'Hx.y)  = u(x)  J sech2(^j)dy  + UQ  f dy  + 5(  -D ) , 

-D  -D 

and 

D D 

5(D)  = U ( x ) J sech2(a-^-y)dy  + UQ  J dy 
-D  -D 


+ 5(-D)  « C2.  (2.22) 
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For  a symmetric  jet  profile  about  y = 0,  C-|  = C2-  Equations 
(2.21)  and  (2.22)  are  conveniently  simplified: 


and 


’Kx.y)  = U ( x ) d ( x ) {tanh(j^j)  + tanh(j^y)} 


+ UQy  + UqD  + J(-D), 


5(D)  = 2U  ( x ) d ( x ) tanh(j^y)  + 2U„D  + •$  ( - D)  . 


(2.23) 


From  (2.23),  the  jet  velocity  at  y=0  is  given  by 


(2.24) 


U(x)  = 


■5(D)  - 2UqD  - *(-D) 


coth  (JTTT1  • 


2ITY) 

Therefore,  if  d(x)  varies  slowly  in  x,  then,  so  does  the 
basic  flow.  The  x-variation  for  the  characteristic  length 
scale  is  gi ven  by 


(2.25) 


( ( x - x n ) / 

i 850  km  + 350  km  <.  cos[2t  — ] > , x ^ x 

‘ I L \ 
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1 200  km 


- I 


• x 1 xo  ) 


. (2.26) 


Here  x is  the  longitude  where  the  x-variation  of  the  cosine 
function  starts  and  L is  the  wavelength  of  this  variation. 

Figures  1-3  show  the  basic  fields  of  streamfuncti on 
5(xiy)»  zonal  velocity  u ( x ,y ) and  vorticity  ;(x,y), 
respectively.  Here  L is  set  to  43,000  km.  The  domain  is 
40,125  km  long  and  4,000  km  wide  (D  = 2,000  km).  At  y=0, 
the  basic  zonal  velocity  is  -15.51  m sec'1  at  the  inflow 
boundary  (x  = 14,625  km)  increasing  slowly  to  a maximum  value 
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of  -30  m sec-1  at  x=0.  From  this  longitude  the  central 
velocity  slowly  decreases  downstream  to  a minimum  value  of 
-13.4  m sec-^  at  x = -21,375  km.  Between  this  longitude  and 
the  outflow  boundary  (x  = -25,500  km),  the  basic  flow  is 
parallel.  From  Eq.  (2.26),  we  note  that  the  ch a r ac te r i s t i c 
length  scale,  d(x),  varies  between  500  km  at  ^(x*°)max  and 
1 200  km  at  u ( x , o ) mi n . 
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Ill w LOCAL  STABILITY  OF  BASIC  FLOW 

The  local  growth  rate  of  the  variable  mean  flow  of  the 
numerical  model  is  first  determined  by  a parallel  flow  (i.e., 
no  x-variation)  numerical  model  of  Williams  e_t  aj_.  (1971)  in 
order  to  gain  some  insight  on  the  stability  characteristics 
of  the  mean  flow.  This  model  is  hereafter  referenced  as  the 
parallel  flow  model.  By  setting  3-j )/3x  = 0 in  (2.9)  of  the 
numerical  model,  the  governing  equation  of  this  model  is 

[ft  + =ll+  * C*  - 4]  If  ■ ° • (3.1) 

°y 

The  symbols  of  (3.1)  have  the  same  meaning  as  the  numerical 
model  of  Chapter  II,  except  that  the  basic  flow  is  given  by 

u(y)  = - U sech2(^-)  , (3.2) 

where  U is  a specified  constant  that  scales  the  magnitude  of 
the  central  velocity  of  the  Bickley  jet  (y=0).  The  character- 
istic length  d is  also  a specified  constant.  Assuming  that 
all  perturbation  quantities  are  periodic  in  x,  Eq . (3.1)  is 
finite  Fourier  transformed  in  x with  wavenumber  k,  and  is 
solved  as  an  initial  value  problem.  This  approach  gives  the 
phase  speed,  growth  rate,  and  the  wave  structure  of  the  most 
unstable  mode.  It  is  convenient  to  write  in  the  following 
form : 

<J>'(x,y,t)  = A (y  , t ) cos(kx)  + B(y,t)  sin(kx).  (3.3) 
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Equation  (3.3)  is  substituted  into  (3.1),  and  the  coefficients 
of  coskx  and  sinkx  are  separately  set  equal  to  zero.  The 
cos kx  term  gives: 

[- 
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and  the  sinkx  term  gives 
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(3.4) 


(3.5) 


(3.6) 


The  boundary  conditions  are 

<|i 1 ( x , D , t ) = ip  ' ( x , - D , t ) = 0. 

Equations  (3.4)  and  (3.5)  are  written  in  finite  difference 
form  such  that  the  second  derivative  with  respect  to  y of  a 
typical  variable  A is  approximated  as 


A . 


C1^]  = 


- 2 A . + A.  , 
J UlL 


] , 


(3.7) 


ay 


j 


Ay 


where  j is  the  y-grid  index  and  Ay  is  the  distance  between 
grid  points.  The  model  has  32  grid  intervals  and  Ay  is  set 
equal  to  125  km  (the  width  of  the  channel  is  4000  km).  This 
y-grid  structure  is  the  same  one  used  in  the  complete 
numerical  model.  Centered  time  differences  are  used  for  all 
quantities  except  those  involving  friction.  The  time  step 


(At)  is  set  to  1 hour.  The  friction  terms  are  evaluated  at 
the  previous  time  step  in  order  to  insure  linear  computa- 
tional stability  (Haltiner,  1971).  The  integration  begins 
with  a forward  time  step.  The  boundary  conditions  (Eq.  3.6) 
become 

A = B = 0 at  y = +D.  (3.8) 

Equations  (3.4)  and  (3.5)  are  solved  for  the  tendencies  by 
the  exact  method  of  Richtmyer  (1967).  These  equations  can  be 
solved  numerically  as  a function  of  3,  D^,  u(y),  and  k for 
any  initial  conditions.  These  equations  are  integrated  until 
their  solution  becomes  exponential  in  behavior.  The  basic 
flow  u(y)  is  specified  for  a number  of  selected  longitudes 
of  the  variable  mean  flow  of  Fig.  2.  The  wavenumber  k is 
also  specified  for  each  integration. 

In  general,  these  equations  have  a set  of  discrete  normal 
mode  solutions  as  well  as  a continuous  spectrum  of  solutions 
(Case,  1960,  Pedlosky,  1964  and  Yanai  and  Nitta,  1968).  Only 
the  normal  mode  solutions  can  give  significant  growth  and  the 
most  unstable  mode  will  dominate  after  a sufficient  period 
of  time. 

This  parallel  flow  model  is  numerically  integrated  to 
150  days.  The  initial  disturbance  amplitude  has  a north- 
south  structure  with  no  tilt.  The  ei gensol uti on  obtained 
includes  wave  structure  of  the  most  unstable  mode,  and  its 
growth  rate  and  phase  speed.  By  selectively  testing  different 
wavelengths,  the  most  unstable  wavelength  is  determined. 


Figures  4-7  illustrate  the  ei gensol uti ons  of  the  most 
unstable  discrete  mode  for  the  selected  jet  profiles.  Since 
the  jet  flow  and  the  wave  disturbance  are  both  symmetric 
about  y = 0,  only  the  lower  half  of  the  y domain  is  shown. 

It  is  interesting  to  note  that  the  eigensolution  amplitude 
has  three  maxima,  one  at  y=0  and  one  on  each  wing  (only  lower 
wing  is  shown)  of  the  jet  approximately  600-800  km  from  y=0. 
It  is  apparent  that  when  the  jet  is  relatively  sharp  (i.e., 
d < 700  km),  the  central  maximum  predominates.  When  the  jet 
is  relatively  smooth  and  broad  (i.e.,  d > 800  km),  the 
maximum  on  the  wings  predominates.  Figure  8 shows  the 
growth  rate  (n)  corresponding  to  the  most  unstable  wave- 
lengths (L)  as  a function  of  x based  on  the  parallel  flow 
model.  In  a local  stability  sense,  the  most  unstable 
wavelengths  range  between  3650  km  at  x=0  where  the  jet  has 
its  maximum  velocity  and  4600  km  near  the  inflow  and  outflow 
regions  of  the  jet  regime.  We  observe  that  locally  the 
greatest  instability  based  on  both  tilt  and  growth  rate  is 
indicated  where  the  jet  achieves  its  maximum  central  velocity 
(u  = -30  m sec'1).  This  parallel  flow  model  can  only  give 
the  most  unstable  discrete  mode;  it  can  not  depict  any 
dynamic  damping.  This  is  reflected  in  Fig.  7(b)  which 
corresponds  to  x = -22,500  km,  outside  the  unstable  region. 
Here  the  solution  exhibits  essentially  no  tilt  and  the  growth 
rate  asymptotically  approaches  the  linear  frictional  damping 
rate.  This  is  further  discussed  in  Chapters  VI  and  VII. 
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Parallel  flow  model  wave  structure  of  the  most  unstable  mode 
for  longitudes:  (a)  x = +7500  km,  and  (b)  x = +9375  km. 
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Parallel  flow  model  wave  structure  of  the  most  unstable  mode 
for  longitudes  (a)  x = +14625  km,  and  (b)  x = -22500  km. 
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These  results  are  consistent  with  those  calculated  by 
Kuo  (1973)  in  a numerical  study  of  a parallel  Bickley  jet 
and,  in  a local  stability  sense,  may  be  viewed  as  a first 
approximation  of  the  behavior  of  moving  waves  within  a 
zonally-varying  mean  wind. 
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IV.  FINITE  DIFFERENCE  EQUATIONS 


The  governing  Eq.  (2.9)  of  the  complete  numerical  model 
is  written  with  the  advective  terms  in  Jacobian  form: 


,23- 


3 t 


- J(p,vV  ) - J(p'  ,v2p) 


~ o y» 

3" 


3 x 


- Df7  S 


(4.1) 


where 


J(a,b)  = 


aa  3b  3a  3b 


(4.2) 


3 x 3y  3y  3 x 

for  any  two  scalar  quantities  a and  b.  The  Arakawa  (1966) 
finite  difference  approximation  for  the  Jacobian  is  used. 

This  scheme  conserves  both  the  mean  square  vorticity  and  the 
mean  kinetic  energy. 

We  define  the  tendency  of  the  disturbance  streamfunction 
by 


T ^ 


i , J 


= rii— ] 

L3t  J 


(4.3) 


^ , J 


where  superscript  t is  the  current  time,  and  subscripts  i 
and  j refer  to  the  x and  y grid  points,  respectively.  In  the 
following  finite  difference  equations,  At  is  the  time  step 
and  is  set  equal  to  one  hour.  Ax  and  Ay  are  the  x and  y 
grid  point  intervals,  respectively.  Equation  (4.1)  is  written 
using  the  Matsuno  or  Eul er-backward  finite  difference  scheme: 
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(4.5) 
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The  finite  difference  form  for  the  Laplacian  of  the  tendency 
is  given  by:  (4.8) 
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and  for  the  Jacobian  terms,  it  is  given  by: 
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where 
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C'/,72^]  is  expanded  in  the  same  manner  as  for 


i . J 


(4.12) 
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The  Poisson  Eq.  (4.1)  is  solved  for  the  tendency  with 
a direct  method  developed  by  Sweet  (1971).  This  direct  method 
solves  a finite  difference  approximation  to  Poisson's  equation 
on  a rectangular  domain  with  Dirichlet  boundary  conditions. 

The  finite  difference  form  of  the  boundary  conditions  are  as 
f o 1 1 ows : 

a)  Northern  (j  = J)  and  southern  ( j = 0 ) boundaries: 


= 0 : 


i =0 , 1 for  j = 0 and  j=J  , 


(4.13) 
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b)  Inflow  or  eastern  boundary  (i=I,  j=0,  j ) 
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where  u>  is  the  specified  frequency  and  k T . is  a specified 
wavenumber  structure.  Time  is  discretized  by: 

t = «.At,  where  a = 0,1, 2, 3, 4 ...  (4. 

A.,  and  B.  are  the  y-structure  coefficients  of  the  forcing 
j , J 

(4.14).  These  coefficients  were  determined  from  the  eigen 

solution  of  the  parallel  flow  model.  The  velocity  profile 

for  determining  A-  and  B.  of  the  inflow  boundary  condition 

J J 

was  selected  from  the  inflow  region  of  the  numerical  model 
This  profile  is  described  in  Chapter  V . 

c)  Outflow  or  western  boundary  (i=0,  j = 0 , J) 


I + A X / ( C 


1 + A X / ( C 


* Vj  «■ 


(4.22) 


The  tendency  Tg  ^ ^ (4.20)  and  s'g^j  (4.21)  were  derived  from 
the  analytical  Sommerfeld  radiation  condition  (2.16)  and 
(2.17),  respecti vely . Since  (2.16)  is  not  in  the  usual  form 
of  the  Dirichlet  boundary  condition,  the  direct  Poisson 
equation  solver,  POISDD,  had  to  be  modified.  Details  of  this 
modification  are  given  in  the  Appendix. 
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V.  EXPERIMENT  I 

The  investigation  of  the  numerical  model  of  the  variable 
mean  wind  comprises  two  experiments.  Experiment  I is  the 
principal  one  and  will  be  discussed  in  detail  in  this  chapter. 
In  Experiment  II  a constant  velocity  is  added  to  the  mean 
flow  and  its  results  will  be  discussed  in  Chapter  VII. 

The  forecast  equation  is  integrated  in  time  from  an 
initial  state  of  ' = 0.  Tne  periodic  forcing  on  the  eastern 
boundary  causes  the  interior  s t re amf un ct i on  to  grow  and  the 
integration  is  continued  until  the  time  variation  is  periodic 
everywhere  with  the  forcing  frequency.  By  H?v  70.  t.hp  model 
has  achieved  the  fully  periodic  state  and  the  wave  packet 
envelope  has  become  quas i -stati onary  . The  forcing  frequency 
is  varied  until  the  value  which  gives  a maximum  perturbation 
ampl i tude  i s found . 

The  following  values  are  used  in  Experiment  I: 

2D  = 4000  km,  x^-x^  = 4 01  25  km,  u(0,0)  = - 30  m sec~\ 

U0  = 0,  L = 43000  km.  Ax  = 375  km,  Ay  = 125  km. 

For  these  parameters  the  maximum  response  occurs  for  a 
forcing  period  of  3.25  days.  It  requires  several  complete 
integrations  to  refine  the  inflow  and  outflow  boundary 
conditions  for  the  specified  forcing  period.  Wavenumber  (k) 
for  the  inflow  vorticity  boundary  condition  (4.13)  was 
determined  by  observing  the  predominant  wavelengths  near  the 
inflow  region.  The  radiation  phase  velocity  (c r ) in  the 
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outflow  boundary  condition  (4.15)  was  similarly  determined 
by  observing  the  predominant  wavelengths  near  the  outflow 
region.  Since  the  domain  is  periodic,  the  phase  velocity  is 
approximated  by  the  simple  relation: 


Here  oj  is  the  local  frequency  which, in  the  fully  periodic 

state, is  equal  to  the  forcing  frequency,  and  kQ  is  the  wave- 

number  measured  near  the  outflow  region.  Thus,  for  Experiment 

I,  a phase  velocity  (c  ) of  -13  m sec"1  is  used  for  the 

' o 

outflow  Sommerfeld  radiation  boundary  condition. 

A.  RESULTS 

The  solution  for  the  Case  I experiment  becomes  fully 
periodic  after  70  days.  The  field  at  t=70  days  is  shown 
in  Fig.  9.  An  entire  train  of  barotropic  waves  actually 
exists  throughout  the  length  of  the  channel,  but  the  waves 
upstream  of  x = -4,500  km  and  in  the  outflow  region  are  not 
shown  in  Fig.  9 because  of  their  relatively  small  amplitude. 
The  maximum  wave  amplitude  occurs  at  x = -12,750  km  and  is 
2 and  4 orders  of  magnitude  larger  than  at  the  jet  maximum 
(x=0)  and  the  inflow  boundary  (x  = 14,625  km),  respectively. 
Figure  10  shows  the  lower  half  of  phase  angle  tilt  s*(x,y) 
of  the  wave  disturbance  y-structure  for  various  values  of 
longitude  x.  We  note  that  the  waves  upstream  of  approximately 
x = -20,000  km  tilt  opposite  to  the  mean  wind  shear,  which  is 
necessary  for  barotropic  instability.  In  fact  the  tilt  near 
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the  inflow  boundary  shows  relatively  weak  instability  and 
this  tilt  slowly  increases,  reaching  a maximum  at  x=0  where 
the  jet  velocity  is  maximum.  Further  downstream,  the  tilt 
slowly  decreases  and  eventually  reverses  near  x = -20,000  km. 
This  behavior  is  consistent  with  the  growth  rates  for  the 
most  unstable  wavelengths  shown  in  Fig.  8 which  were  computed 
using  the  parallel  flow  model.  Near  the  outflow  boundary, 
however,  the  tilt  is  reversed.  This  indicates  dynamic 
stability  or  a flow  of  energy  from  the  disturbance  back  to 
the  mean  flow.  Therefore,  we  observe  that  the  tilt  of  the 
wave  disturbance  qualitatively  adjusts  to  the  local  stability 
of  the  mean  flow  (see  Figs.  4-7).  We  recall  that  the  parallel 
flow  model  can  only  solve  for  the  most  unstable  discrete  mode, 
thus,  dynamic  damping  is  not  indicated  in  Fig.  7(b)  for  the 
outflow  region. 

Figure  11  shows  the  envelope  of  the  wave  packet, 
evaluated  at  y = -750  km,  where  the  disturbance  amplitude  is 
large.  This  envelope  <-j'(x)>  is  obtained  by  recording  the 
maximum  and  minimum  y ' values  that  occur  at  each  longitude 
over  a period  of  ten  days  after  the  solution  becomes  fully 
periodic.  If  a larger  time  interval  were  chosen,  the  envelope 
would  not  change.  Note  that  the  maximum  amplitude  occurs  in 
the  area  where  the  local  growth  rate  becomes  zero  (see  Fig.  8). 
The  smoothness  of  the  fields  in  Figs.  9 and  11  indicate  that 
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The  iji 1 fields  for  days  83  through  86  are  shown  in  Figs. 

12  through  15,  respectively,  to  illustrate  a time  sequence 
of  the  -y  barotropic  wave  train.  A periodicity  of  3.25  days 
can  be  determined  by  examining  the  time  series  at  each  grid 
point  (not  shown).  We  also  note  from  this  sequence  that  there 
are  three  maxima  of  <p ' in  the  wave  structure.  One  maximum  is 
at  y=0  and  another  is  on  each  wing  of  the  jet  near  y = +_7 5 0 km. 
The  maximum  on  the  wing  clearly  predominates.  Figure  16  shows 
the  variation  of  wavelength,  L(x),  for  y=0  and  y = +750  km 
corresponding  to  the  latitudes  of  the  observed  three  maxima 
of  t|/ 1 . The  disturbance  wavelength  is  about  4100  km  initially 
near  the  inflow  bO|,nd»^v  Thp  wavplpnnfh  npar  v = +750  km 
is  larger  upstream  and  smaller  downstream  of  x=0  than  at  y=0. 

At  y = +750  km,  the  maximum  wavelength  of  5060  km  occurs 
ap proxi ma te 1 y 950  km  upstream  of  x=0  while  the  minimum  wave- 
length of  3600  km  occurs  near  the  outflow  boundary.  At  y=0, 
the  maximum  wavelength  of  4900  km  occurs  about  950  km  down- 
stream of  x=0  while  the  minimum  wavelength  of  about  4000  km 
occurs  near  x = -13,500  km  with  the  wavelength  increasing 
slowly  further  downstream  to  the  outflow  boundary.  The 
range  of  the  wavelength  is  between  approximately  3600  km  and 
5000  km,  which  is  nearly  the  same  range  of  values  obtained 
for  the  most  unstable  wavelengths  using  the  parallel  flow 
model  in  Chapter  III.  Figure  17  shows  the  disturbance  phase 
speeds,  c *(x),  for  latitudes  y=0  and  y = +750  km  where 
cr*(x)  is  obtained  from 
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Wavelengths  L(x),  from  the  complete  numerical 
udes  (a)  y-0,  and  (b)  y = -750  km. 


Since  the  entire  domain  is  periodic  with  the  forcing 
frequency,  c^*(x)  has  the  same  basic  behavior  as  L(x) 

[Fig.  16]  for  the  corresponding  latitudes.  Namely,  upstream 
of  x=0  the  disturbance  phase  velocity  is  faster  on  the  wings 
of  the  jet  than  at  y=0,  and  downstream  of  x=0  it  is  slower. 
This  behavior  of  both  L(x)  and  cr*(x)  illustrates  the  tilt 
behavior  of  the  waves  as  they  move  into  and  out  of  the  baro- 
tropically  unstable  region  of  the  variable  mean  wind.  More 
results  are  presented  in  the  next  chapter. 

The  parallel  flow  model  is  also  employed  to  determine 
a reasonably  good  forcing  function  that  would  generate 
periodic  wave  disturbances  from  the  eastern  boundary.  This 
is  accomplished  as  follows: 

Equation  (2.14)  defines  the  analytical  form  of  the 
periodic  forcing  which  is  applied  as  an  inflow  boundary 
condition  in  the  complete  numerical  model,  i.e., 


( 0 ,y  , t ) = A (y  ) u>cos  ( wt ) - B(y)asin(at) 


Equation  (4.17)  describes  the  finite  difference  approximation 
to  this  forcing.  The  coefficients  A(y)  and  B(y)  of  (2.14) 
are  obtained  from  the  ei gensol uti on  determined  using  the 
parallel  flow  model  for  the  following  Bickley  jet  profile: 


u (y=0 ) 

d 

L 
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This  profile  is  representative  of  the  mean  flow  near  the 
inflow  region.  Figure  18  shows  the  y-structure  of  these 
coeff i ci ents . 

Two  other  forcing  functions,  which  are  independently 
constructed  from  the  parallel  flow  model,  are  also  tested 
These  functi ons  are  : 


fy-  (0,y,t)  = - a)  cos(— sin(jt),  (5.3) 


and 


3 1 
3 t 


(0,y,t) 


(cos  ( 


s i n ( u>  t ) , 


(5.4) 


where  D is  the  half  width  of  the  channel  (D  = 2000  km). 

The  behavior  of  the  wave  disturbances  for  all  three 
forcing  functions  are  essentially  the  same  except  in  the 
vicinity  of  the  inflow  boundary.  The  parallel  flow  solutions 
for  the  most  unstable  wavelength  show  a weak  decay  rate  in 
the  inflow  region  (see  Fig.  8).  This  is  because  the  linear 
friction  decay  rate  is  slightly  larger  than  the  weak  local 
growth  rate.  We  found  that  this  weak  decay  behavior  near 
the  inflow  boundary  is  achieved  in  the  numerical  model  using 


forcing  function  (2.14).  Use  of  the  other  two  forcing 
functions,  (5.3)  and  (5.4),  realizes  a weak  seconda ry  < * '>  x 
near  the  inflow  boundary.  Therefore,  we  found  through  use 
of  the  parallel  flow  model  that  the  first  forcing  function 
provided  the  better  inflow  boundary  condition  than  the  other 
two,  although  no  substantial  difference  occurs  even  if  the 
other  forcing  functions  are  used. 
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B.  ENERGETICS 

To  understand  the  dynamics  of  barotropic  instability, 
it  is  helpful  to  derive  the  disturbance  kinetic  energy 
equation.  A finite  difference  approximation  form  of  this 
equation  is  used  to  compute  the  energy  budget  of  the  waves. 

The  linearized  momentum  equations  are: 
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After  multiplying  (5.5)  by  u'  and  (5.6)  by  v'  and  adding  the 
resultant  equations,  we  obtain  the  disturbance  kinetic  energy 
equation: 
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(5.7) 


where 


(5.8) 


The  disturbance  fields  are  considered  to  have  a time 
dependence  of  the  form  sin[wt  + e(x,y)].  Thus,  when  (5.7) 
is  averaged  over  one  period  and  over  one  latitudinal  domain 
(-D^y<_0),  the  following  energy  balance  equation  is  obtained: 
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Term  1 is  the  advection  of  k'  by  the  mean  flow  ( u ) . Terms 
2,  3,  4 and  5 are  the  barotropic  exchange  terms  with  the 
mean  flow.  Term  6 represents  the  work  done  by  the  pressure 
field.  Term  7 is  the  dissipation  term  and  represents  an 
energy  sink.  With  the  exception  of  term  6,  all  the  terms  i r\ 
(5.9)  are  determinable  from  the  J and  ' fields.  Thus,  the 
pressure  work  term  is  determined  as  a residual  in  Eq.  (5.9). 
The  energy  equation  is  evaluated  after  the  complete  numerical 
model  achieves  the  fully  periodic  state. 

Figure  19  shows  the  energy  balance  as  function  of 
longitude.  Clearly,  the  maximum  value  of  each  of  the  terms 
in  the  energy  equation  occurs  downstream  of  x=0.  Curve  (a) 
represents  the  Reynold's  stress  term,  <( - u ' v ' 3 u/ 3y ) > , where 
the  maximum  value  occurs  approximately  at  x = -10,975  km. 
Curve  (b)  represents  the  advection  term,  <(-u3k'/3x)>  , of 


k uy  une  l u 1 1 a i current.  Curve  (c)  represents  the  linear 

dissipation  term  <(2D^k')>.  We  observe  that  the  largest 

dissipation  takes  place  near  the  maximum  amplitude  of  the 

wave  packet  envelope,  <^'(x)>  . Curve  (d)  represents  the 

ma  x 


pressure  work  term,  <C  ^)>  It  also  has 


a maximum 


near  < ^ 1 ( x ) > 


max 


The  remaining  barotropic  terms  are  not 
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turb 


included  in  Fig.  18  because  they  are  2 to  3 orders  of 
magnitude  smaller  than  <(-  u'v'3u/3y)>. 

Recall  that  for  Experiment  I,  <^'(x)>  occurs  approxi- 

max 

mately  at  x = -12,750  km.  We  observe  that  the  advection  of 

k1  by  u plays  a dominant  role  in  the  energy  balance.  Once 

the  fully  periodic  state  is  achieved,  the  generation  of 

disturbance  energy  by  barotropic  instability  essentially 

occurs  in  the  region  between  the  <-p  1 > (x  = -12750  km)  and 

max 

the  jet  maximum  (x=0).  It  is  in  this  region  that  the 
strongest  positive  correlations  occur  between  u ' v ' and 
- 3 u / 3 y in  the  Reynold's  term.  This  is  consistent  with  the 
observed  behavior  of  the  waves  (see  Figs.  10  and  11)  since 
this  strong  correlation  of  u ' v 1 is  dependent  on  maximum 
amplitude  of  ip  1 squarred  and  the  maximum  tilt  of  the  wave 
structure  with  the  shear  of  the  mean  flow  -3u/3y.  It  is 
also  from  this  source  region  that  u transports  a large  part 
of  this  generated  k'  downstream  to  a point  centered  near 
x = -1  5 ,000  km. 

Strongest  dissipation  of  k ' occurs  in  the  region  centered 
near  < ip ' ( x ) >ma  v . It  is  in  this  region  that  the  strongest 
spatial  gradients  of  tp  1 occur.  The  largest  values  of  pressure 
work  also  occur  in  this  region  where  the  pressure  work  term 
is  computed  as  a residual.  We  note  that  its  location  of 
occurrence  is  consistent  with  the  large  amplitudes  and  spatial 
gradients  of  ip'  of  the  region.  We  also  note  that  <(-u3k'/3x)> 


L 


63 


vanishes  at  <4i'(x)>max  since  k1  is  maximum  at  this  longitude. 
The  disturbance  energy  balance  which  is  shown  in  Fig.  19 
graphically  illustrates  the  basic  energetics  involved  in 
supporting  the  large  amplitudes  and  spatial  gradients  of  -h 1 
a considerable  distance  downstream  of  the  jet  maximum. 

Clearly,  the  Reynold's  stress  term,  < ( u ' v 1 3u/3y )> , is  the  only 
dominant  source  for  disturbance  kinetic  energy. 
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VI.  SIMPLE  MECHANISTIC  ANALYTICAL  MODEL 


A.  FORMULATION 

In  this  chapter  we  develop  a simple  mechanistic  analyti- 
cal model  which  uses  the  locally  determined  parallel  flow 
solutions.  The  results  of  this  simple  model  can  be  compared 
with  those  of  the  complete  numerical  model  to  obtain  some 
measure  of  the  adjustment  of  the  transient  barotropic  waves 
to  the  local  stability  of  the  variable  mean  wind. 

The  following  equation  allows  for  propagation  and  growth 
or  decay: 

fr* + cr(x)  rf  = n(*h'  * (6J) 

where  -y  represents  the  disturbance  s treamf uncti on . c (x) 
is  the  local  phase  velocity  and  n(x)  is  the  local  growth  rate. 
We  observe  that  if  cr  and  n are  independent  of  x,  then 
Eq.  (6.1)  is  exact,  and,  in  general,  should  give  a reasonable 
approximation  to  the  downstream  variation  of  ' . 

Let  us  consider  the  case  where  ^ ' has  a solution  of  the 

form : 

*'  = F'(x)eiwt,  (6.2) 

where  u i s a specified  frequency.  After  substituting  (6.2) 
into  (6.1),  we  obtain 

1uF'  + cr(x)  = n ( x ) F ' . (6.3) 
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Rearranging  (6.3),  we  get 


1 dT_  n ( x ) - i oi 

F 1 dx  cr(x)  ’ 

which  has  a solution  of  the  form  of 


x 


F ' (x) 


p;  exp 


-<  / 


0) 


dx>  exp 


(6.4) 


(6.5) 


.( 

x„ 


(6.6) 


n ( x ) 
c^Jx) 


Some  fundamental  insight  can  be  gleaned  from  (6.6).  Here 
the  amplitude  of  p ' must  be  specified  at  some  initial  point, 
xQ,  which  could  be  the  inflow  point.  The  local  wavenumber 
is  «/cr  and  the  spatial  growth  rate  is  n/c^.  These  are  the 
real  and  imaginary  parts  of  the  wavenumber,  if  we  were  to 
write  ip ' = A exp  (ikx-iwt). 

In  this  simple  analytical  model,  we  are  interested  in 
the  growth  and  decay  of  the  wave  packet.  If  the  oscillatory 
behavior  of  Eq.  (6.6)  is  dropped,  the  following  equation 
expresses  the  envelope  of  the  wave  packet: 


i-  (x) 


i-  exp  y 
0 I 


I 


n ( x ) 
cr(xT 


dx 


\ 


(6.7 


Equation  (6.7)  gives  the  exponential  growth/decay  behavior 
of  a wave  disturbance  as  it  travels  within  its  wave  packet 
envelope.  This  behavior  is  determined  by  the  integral 
effects  of  the  local  stability  properties  of  the  mean  flow. 
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B.  APPLICATIONS 

Figure  20  compares  the  envelope  for  1 ( x ) from  the 
complete  numerical  solution  (Fig.  11)  with  the  F(x)  given 
by  Eq.  (6.7)  for  two  values  of  xQ.  In  this  figure  the 
amplitude  is  plotted  against  upstream  distance  from  the 
initial  point  xq.  The  analytical  model  uses  the  locally 
determined  parallel  flow  growth  rates  (n)  and  Dhase  speeds 
Ccr)  based  on  the  wavelengths  measured  from  the  complete 
numerical  solutions.  Therefore,  Eq.  (6.7)  is  applied  at 
latitude  y = + 750  km  where  the  wavelengths  are  measured. 
The  curve  for  F(x)  can  he  shifted  up  or  down  by  changing 
xq,  but  its  shape  does  not  change.  In  Figure  20  all  the 
curves  have  a maximum  amplitude  at  x=  - 12,750  km  where 
the  local  growth  rate  is  zero;  however,  the  simple  model 
has  a lower  maximum  growth  rate  and  slower  damping  rates 
than  the  complete  numerical  model. 

The  lower  portion  of  Fig.  21  contains  the  wavelength 
measured  from  the  numerical  solution  at  y = - 750  km.  In 
the  upper  portion  of  the  diagram  are  the  phase  velocities, 
c„  from  the  parallel  flow  model  and  c * from  the  complete 
numerical  model  for  the  same  latitude.  Here  c^,  is  computed 
with  the  use  of  the  wavelength  shown.  The  two  phase  speeds 
have  similar  behavior  although  c^,*  is  shifted  slightly 
upstream  and  has  larger  variations.  We  also  note  that  c^* 
is  generally  smaller  than  cr  except  for  the  region 
-1 500  km  < x < 6750  km. 
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Similarly,  Fig.  22  contains  the  wavelengths  and  phase 

speeds  for  latitude  y = 0.  Again,  the  two  phase  speeds  have 

similar  behavior,  but  at  this  latitude  c * is  shifted 

r 

slightly  downstream  and  the  differences  between  phase  speeds 
is  generally  less  than  at  y = -750  km.  Again  we  note  that 
c^*  is  generally  smaller  than  c.  except  for  the  region: 

-2250  km  < x < 2250  km.  This  indicates  that  both  the  wave- 
length L(x)  and  the  phase  velocity  cr*(x)  between  latitudes 
y=0  and  y = -750  km  are  maximum  where  the  mean  flow  is 
s tronges  t . 

Figure  23  contains  tne  spatial  growth  rates,  m and  m* , 
from  the  parallel  flow  model  and  from  the  complete  numerical 
model,  respecti vely . For  the  parallel  flow  model,  we  obtain 

m = - n/cr  . (6.8) 

The  value  of  m*  is  computed  directly  from  the  envelope.  We 
note  in  Fig.  23  that  m*  has  a larger  maximum  than  m and  the 
maximum  is  shifted  slightly  downstream  from  the  jet  maximum. 
Both  curves  pass  through  0 at  x = -12,750  km  which  is  the 
maximum  for  both  wave  packet  envelopes.  However,  m* 
shows  a much  larger  damping  in  the  outflow  region  which  can 
also  be  seen  in  Fig.  20.  In  fact  the  parallel  flow  solution 
damps  at  the  rate  given  by  the  frictional  coefficient 
divided  by  the  phase  speed.  This  is  expected  in  the  baro- 
tropically  stable  outflow  region  when  using  the  parallel 
flow  model  because  this  model  can  only  give  the  eigen- 
solution  of  the  most  unstable  discrete  mode. 
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I TUDE  : 


22.  Experiment  I.  Latitude  y = 0.  (a)  phase  velocities, 

cr(x)  from  parallel  flow  model,  and  cr*(x)  front  comple 


rates,  m(x)  from  the  parallel 
complete  numerical  model. 


The  solutions  observed  in  the  numerical  model  (see 
Fig.  10)  toward  the  outflow  boundary  tilt  in  the  same  sense 
as  the  mean  wind  shear  which  gives  dynamic  damping.  This 
tilt  behavior  is  also  indicated  by  the  divergence  of  phase 
speed  in  the  outflow  region  between  y=0  and  y = -750  km 
(see  Figs.  21  and  22).  This  dynamic  damping  appears  to 
be  due  to  continuous  spectrum  solutions  (Case,  1960;  Yanai 
and  Nitta,  1968)  which  are  not  included  in  the  normal  mode 
solutions  that  are  employed  in  the  simple  integral 
[Eq.  (6.6)].  In  fact  we  observe  that,  in  the  outflow  region, 
the  dynamic  damping  is  1/t  or  faster,  increasing  near  the 
outflow  boundary.  In  Fig.  23  the  m*  curve  is  skewed  slightly 
to  the  left  with  respect  to  m and  the  jet  maximum.  This 
could  be  the  result  of  the  tilt  structure  in  the  wave  lagging 
spatially  behind  the  expected  value  from  the  local  stability 
conditions.  This  would  give  a smaller  growth  rate  on  the 
upwind  side  of  the  jet  maximum  and  a larger  growth  rate  on 
the  downwind  side.  This  effect  is  indicated  by  the  asymmetry 
of  the  m*  curve  with  respect  to  the  m curve  in  the  unstable 
region  of  the  mean  wind.  However,  this  is  not  the  only  effect 
involved  because  the  most  striking  feature  of  Fig.  23  is 
the  fact  that  the  maximum  value  of  m*  is  significantly 
larger  than  the  maximum  value  of  m.  We  have  already  shown 
that  the  most  important  disturbance  energy  production  term 
(see  Fig.  18)  in  the  energy  computations  with  the  numerical 
solutions  is  proportional  to 


i. 
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< ( -u ' v ' 3u/ 3y ) > . 


This  term  is  the  only  source  term  for  the  parallel  flow 
model,  and  it  depends  on  the  phase  tilt  in  the  disturbance 
field.  Figure  24  compares  the  phase  tilt  of  the  numerical 
model  solution,  9*,  and  that  of  the  most  unstable  wavelength 
of  the  parallel  flow  model,  9 , at  x = 0,  x = 3750  km,  and 
x = 13125  km.  It  is  evident  that  at  the  jet  maximum  the 
tilt  from  the  numerical  model  is  significantly  larger  than 
for  the  tilt  from  the  parallel  flow  model.  This  is  also 
evident  but  to  a lesser  degree  at  x = 3750  km.  These  results 
are  consistent  with  the  larger  growth  rates  for  the  variable 
jet  flow.  This  stronger  tilt  is  also  indicated  by  the 
comparison  of  the  c^*  and  c^  curves  in  the  region  near  and 
slightly  upstream  of  the  jet  maximum  for  latitudes  y=0  and 
y = -750  km  (see  Figs.  21  and  22).  We  note  that  at 
y = -750  km,  c^,*  is  relatively  larger  than  cr  near  the  jet 
maximum.  This  difference  is  significantly  reduced  and 
shifted  slightly  downstream  at  y=0.  We  also  note  that  near 
the  inflow  region  (x  = 13125  km),  the  phase  tilt  of  the 
numerical  solution  lags  slightly  behind  the  oarallel  flow 
solution  (see  Fig.  24c).  Apparently  the  downstream  varia- 
tion of  u augments  the  phase  tilt  which  gives  a larger 
growth  rate. 

In  comparing  the  parallel  flow  theory  and  the  complete 
numerical  model,  a source  of  error  could  be  introduced 


through  the  zonal  resolution.  This  is  because  in  the  former 
the  zonal  variation  is  specified  by  a SDectral  representa- 
tion, while  in  the  latter  it  is  resolved  by  finite  differences. 
Since  the  wavelength  appearing  in  the  complete  numerical  model 
varies  between  10  and  13.5  ax,  such  resolution  normally  gives 
a good  approximati on  to  the  exact  solution  (Haltiner,  1971). 
Thus,  this  error  should  not  be  very  significant.  In  order  to 
confirm  this  we  carried  out  an  experiment  in  which  a parallel 
mean  flow  is  specified  by  the  full  two-dimensional  model 
finite  differences  of  the  complete  numerical  model.  The 
growth  rate  and  phase  speed  calculated  from  the  experiment 
agree  very  well  with  those  obtained  by  the  semi - s pe ct r a 1 
parallel  flow  model,  the  difference  being  only  1-2".  Thus 
we  may  conclude  that  the  differences  between  the  numerical 
solution  and  the  local  solution  of  the  parallel  flow  theory 
are  genuine  and  are  not  due  to  the  differences  in  the  model 
resolution. 

C.  EFFECTS  OF  LINEAR  FRICTION 

Figure  23  shows  that  the  local  growth  rate,  determined 
by  the  parallel  flow  model,  changes  sign  at  x = -12,750  km 
and  that  it  approaches  the  linear  friction  decay  rate  near 
the  outflow  region.  We  also  note  that  the  maximum  aplitude 
of  both  the  analytical  and  numerical  model  wave  packet 
envelopes  occurs  at  the  lontitude  where  n(x)  vanishes. 

Based  on  these  results,  <^‘>  is  expected  not  only  to 
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decrease  in  magnitude  but  to  occur  further  upstream  if  the 


linear  frictional  rate  is  increased.  Likewise,  the  opposite 

effect  would  be  expected  if  this  friction  were  decreased. 

This  hypothesis  is  tested  by  re  running  Experiment  I using 

- 5 - 1 

a larger  friction  of  = .25x10  sec  which  is  equiva- 
lent to  e-folding  decay  time  of  4.63  days.  We  found  that 

<\p‘>  , # occurred  approximately  1 700  km  further  upstream  and 
max 

with  an  amplitude  of  an  order  of  magnitude  smaller  (.4  vice 
5.9)  than  for  the  original  experiment. 
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VII.  EXPERIMENT  II 


In  Experiment  II  a constant  advective  speed  of 
Uo  = -5  m sec'^  is  added  to  the  mean  flow,  otherwise  the 
parameters  are  the  same  as  those  used  in  Experiment  I.  We 
hope  this  modification  will  shed  some  light  on  the  adjustment 
process  of  a faster  moving  wave  chain.  As  in  Experiment  I, 
the  analytical  model  also  uses  locally  determined  parallel 
flow  growth  rate  (n)  and  the  phase  speed  (c)  based  on  the 
wavelength  measured  from  the  complete  numerical  solution. 

The  local  phase  speed  c is  now  related  to  the  Doppler- 
shifted  phase  speed  c by  c = c + 1)  . The  maximum  -p ' response 
in  this  experiment  is  obtained  with  a forcing  period  of  2.5 
days.  Since  the  Dopp 1 e r- s h i f te d frequency  may  be  written 
as  a = k(c*  - U ) = u - kllo  , this  shorter  period  actually 
gives  approximately  the  same  Doppler-shifted  frequency  as 
that  of  Experiment  I (where  the  period  is  3.25  days  with 
UQ  = 0).  Figure  25  shows  the  p ' field  for  day  70. 

Based  on  the  results  of  the  simple  analytical  model 
[Eq.  (6.8)],  we  expect  that  the  constant  advective  velocity 
U will  significantly  reduce  the  local  spatial  growth  rate 

m = - — 1 - (7.1) 


Figure  26  shows  the  spatial  growth  rates  from  the  parallel 
flow  model  (m)  and  from  the  complete  numerical  model  (m*). 


Comparing  this  figure  with  Fig.  23,  the  reduction  in  both  m 
and  m*  from  those  of  Experiment  I is  clear.  This  reduction 
is  evidently  due  to  the  advective  velocity  which  causes  the 
disturbances  of  Experiment  II  to  move  approximately 
-5  m sec’^  faster  than  in  Experiment  I.  Thus,  the  waves 
have  less  time  to  grow  or  damp  in  response  to  the  local 
stability  per  unit  distance  traveled.  The  smaller  spatial 
growth  rate  obviously  affects  the  amplitude  and  spatial 
gradients  of  ip ' . Figure  27  shows  the  wave  packet  envelope 
< 1 ( x ) > of  the  complete  numerical  solution  for  Experiment  II. 
Figure  28  shows  this  envelope  along  with  the  analytical  model 
envelope  F(x)  using  either  xQ  = -375  km  or  xQ  = 3750  km  as 
the  initial  point  (no  discernible  difference  is  noted  between 
the  use  of  these  initial  points).  Both  envelopes  are,  as 
expected,  significantly  smaller  than  those  of  Experiment  I. 

In  addition,  there  is  less  variation  in  the  wavelength 
(3750  km  < x < 4875  km)  compared  to  Experiment  I.  This  is 
also  a result  of  the  faster  phase  speed  which  causes  the  waves 
to  have  less  time  to  adjust  to  the  local  u(x,y). 

Although  the  faster  phase  speed  of  the  waves  in 
Experiment  II  reduces  the  magnitude  of  growth  and  decay  rates 
and  therefore  brings  the  spatial  growth  rate  curves  of  the 
numerical  model  (m*)  closer  to  that  of  the  parallel  flow 
model  (m),  we  notice  a stronger  asymmetry  in  the  m*  curve 
with  respect  to  x=0  in  Fig.  26  compared  to  Experiment  I 
(see  Fig.  23).  In  addition  to  the  maximum  m*  occurring 
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slightly  further  downstream  of  the  maximum  of  m(x=0)  than  in 
Experiment  I,  the  cross-over  points  of  the  m*  and  m curves 
exhibit  more  asymmetry  in  Experiment  II  than  I.  The  first 
cross-over  point  occurs  approximately  at  longitude  x = 8500  km 
and  x = 4000  km  for  Experiments  I and  II,  respectively,  while 
the  second  cross-over  point  occurs  at  approximately 
x = -12000  km  where  the  spatial  growth  rates  vanish.  As  in 
Experiment  I,  we  interpret  this  asymmetry  as  being  due  to  the 
tilt  structure  of  the  waves  lagging  behind  the  expected  value 
from  the  local  stability  condition.  The  fact  that  this  lag 
effect  is  enhanced  by  the  faster  phase  speed  of  the  waves  in 
Experiment  II  can  also  be  seen  in  Fig.  29  where  the  north- 
south  phase  tilt  9*(x,y)  of  the  wave  disturbances  for  both 
experiments  are  compared  at  several  longitudes.  The  unstable 
tilt  in  Experiment  II  is  weaker  upstream  and  stronger 
downstream  of  approxi mately  x = -2000  km  reflecting  the 
stronger  lag  effect  of  the  faster  moving  waves. 

figure  30  compares  the  phase  tilt  of  the  numerical  model 
solution,  9*,  and  that  of  the  most  unstable  wavelength  of 
the  parallel  flow  model,  e , for  two  longitudes  (x=0  and 
x = 3750  km).  Comparing  this  figure  with  Fig.  24  of 
Experiment  I,  it  is  obvious  that  the  apparent  augmentation 
of  the  wave  tilt  due  to  the  zonal  variation  of  u(x,.y)  is  now 
significantly  reduced.  This  behavior  may  at  least  be 
partially  explained  by  the  lag  effect  caused  by  the  advective 
current  U , which  generally  reduces  the  growth  rate  in  the 
vicinity  of  the  jet  maximum. 
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The  lower  portion  of  Fig.  31  shows  the  wavelength 
measured  from  the  numerical  solution  at  y = -750  km.  In 
the  upper  portion  of  the  diagram  are  the  phase  velocities  at 
y = -750  km,  c*  from  the  complete  numerical  model  and  c from 
the  parallel  flow  model  based  on  the  correspond!'  ng  numerical 
solution  wavelength.  The  two  phase  speeds  agree  very  closely 
upstream  of  the  jet  maximum  while  c*  is  smaller  than  c down- 
stream of  this  longitude.  Similarly,  Fig.  32  contains  the 
wavelength  and  phase  speeds  for  latitude  y=0.  At  this 
latitude,  the  two  phase  speeds  have  similar  behavior  through- 
out the  length  of  domain  with  c*  being  app roxi ma t el y 
1 m sec""'  slower  than  c.  Figure  33  shows  the  wavelength 
curves  [L(x)]  for  latitudes  y=0  and  y = -750  km.  At  both 
latitudes  the  maximum  wavelength  occurs  near  the  jet  maximum, 
but  the  wavelength  (and  therefore  phase  velocity)  at 
y = -750  km  is  larger  upstream  and  smaller  downstream  of 
app rox i ma te ly  x = -5000  km.  As  in  Experiment  I,  these 
differences  in  c*(x)  and  L(x)  between  the  two  latitudes  are 
consistent  with  the  observed  tilt  behavior  of  the  waves. 

Other  than  the  foregoing  mentioned  differences,  the 

basic  behavior  of  the  waves  in  Experiment  II  is  similar  to 

that  of  Experiment  I.  This  includes  the  tilt  of  the  waves 

and  its  reversal  from  that  of  growth  to  that  of  dampinq 

downstream  of  <V>  , the  occurrence  of  <-p’>  a consider- 

■ max’  y max 

able  distance  downstream  of  the  jet  maximum,  the  general 
shape  of  the  wave  packet  envelooe,  and  the  longitudinal 
variation  of  the  energy  balance  (Fig.  34). 


X ( 1 000  KM)  - 

energy  balance:  (a)  Reynold's  stress  term, 
ransport  term  <( -uak ' /3x)>;  (c)  dissipation 

1 t) 

ressure  work  term,  <( - - — (u'p')>. 


VIII.  SUMMARY  AND  CONCLUSIONS 

In  this  work  we  studied  the  behavior  of  waves  moving 
into  and  out  of  a barotropi cal ly  unstable  mean  wind  field 
which  varies  both  in  x and  y.  This  mean  wind  field  roughly 
simulates  the  200  mb  easterly  jet  south  of  the  Tibetan  high 
during  the  Northern  Hemisphere  summer.  A rectangular  domain 
is  used  with  a time  periodic  forcing  on  the  inflow  (east) 
boundary  and  a Sommerfeld  radiation  condition  on  the  outflow 
(west)  boundary.  This  allows  the  simulation  of  the  propaga- 
tion of  small  amplitude  waves  through  the  easterly  jet  region. 

The  V 0 r 1 1 C i t V O n M a f i q n n n a ho  f a nlano  i c;  cnl\/oH  with  tho  IK  P 

of  finite  differences,  and,  when  the  boundary  conditions  are 
properly  adjusted,  the  barotropic  waves  move  smoothly  across 
the  region  and  out  the  western  boundary.  After  a certain  time 
period  of  numerical  integration,  the  solution  becomes 
periodic  everywhere  with  the  forcing  frequency  which  is 
specified  on  the  eastern  boundary.  As  the  waves  move 
through  the  jet  regime,  they  grow  or  decay  spatially  in 
response  to  the  local  stability  of  the  mean  flow,  even 
though  at  each  point  the  variation  is  purely  harmonic.  The 
numerical  results  are  compared  with  a simple  mechanistic 
analytical  model  which  is  developed  using  the  local  growth 
rate  and  phase  velocity  computed  from  the  parallel  flow 
theory . 
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Two  experiments  are  carried  out  in  this  investigation. 
The  parameters  used  in  these  experiments  are  the  same, 
except  that  in  Experiment  II  a constant  advective  velocity 
of  UQ  = - 5 m sec'^  is  added  to  the  mean  flow  in  order  to 
gain  a better  understanding  of  the  adjustment  process  of  a 
faster  moving  wave  train. 

The  basic  behavior  of  the  waves  in  both  experiments  is 
similar.  This  includes  the  following: 

1)  The  maximum  amplitude  of  the  waves  occurs  a 
considerable  distance  downstream  from  the  most  unstable  Dart 
of  the  jet  near  where  the  local  growth  rate  of  the  parallel 
flow  theory  vanishes. 

2)  The  tilt  structure  of  the  waves  and  its  reversal 
from  that  of  growth  to  that  of  damping  occurs  downstream  of 
the  maximum  of  the  wave  packet  envelope  in  response  to  the 
local  stability  of  the  mean  flow. 

3)  The  wavelength  (and  therefore  phase  velocity) 
increases  slowly  from  the  inflow  boundary  to  a maximum  value 
near  the  jet  maximum  and  then  decreases  slowly  downstream 
reaching  a minimum  value  near  the  outflow  region. 

4)  The  Reynold's  stress  term  < ( u ' v ' su/ ay  ) > is  the 
dominant  source  for  disturbance  kinetic  energy  and  the 
advective  term  <(-uak!/ax)>  plays  a significant  rule  in 
transporting  a considerable  part  of  this  energy  downstream 
from  the  source  region. 

5)  Large  spatial  damping  occurs  in  the  outflow  region 
due  to  the  presence  of  the  continuous  spectrum  modes. 


93 


r 


In  Experiment  I where  the  mean  wind  vanishes  toward 
large  y,  the  maximum  spatial  growth  rate  in  the  numerical 
model  is  significantly  larger  than  the  spatial  growth  rate 
of  the  parallel  flow  model.  Apparently  the  downstream 
variation  of  the  mean  wind  augments  the  phase  tilt  of  the 
waves  thus  causing  an  increased  growth  rate  as  compared  with 
waves  in  a parallel  flow.  In  Experiment  II  where  the  waves 
move  westward  faster  due  to  the  presence  of  a constant 
- 5 m sec”^  advective  zonal  wind,  the  spatial  growth  rate 
computed  by  the  numerical  model  is  brought  closer  to  that 
given  by  the  parallel  flow  theory,  but  the  asymmetry  of  the 
numerical  growth  rate  curve  with  respect  to  the  jet  maximum 
is  enhanced.  The  reduction  of  spatial  growth  and  decay  rates 
is  mostly  due  to  the  faster  phase  speed  which  causes  the 
waves  to  have  less  time  to  grow  or  damp  per  unit  distance 
traveled.  The  enhanced  asymmetry  is  apparently  due  to  a lag 
effect  wherein  the  adjustment  of  the  waves  is  delayed 
spatially.  The  latter  effect  is  also  reflected  in  the  north- 
south  phase  structure  of  the  waves  and  in  the  smaller  varia- 
tion of  the  wavelength  as  the  waves  propagate  downstream. 

The  augmentation  of  the  phase  tilt  by  the  downstream 
variation  is  also  significantly  reduced  which  may  at  least 
be  partially  explained  by  the  lag  effect. 

Recently  Pedlosky  (1976)  used  a simple  analytical  model 
to  study  the  effect  of  the  downstream  variation  of  a baro- 
clinic  current.  The  mean  flow  in  his  model  changes  abruptly 
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from  a weakly  unstable  regime  with  a constant  growth  rate 
to  a stable  regime  downstream.  His  model  includes  nonlinear 
effects  but  it  does  not  produce  strong  dynamic  damping  due 
to  the  continuous  spectrum  modes.  This  is  perhaps  because 
he  did  not  include  the  lateral  shear  of  the  mean  flow.  Thus, 
Pedlosky  found  that  his  disturbances  remain  substantially 
undiminished  when  they  propagate  from  a locally  unstable 
region  to  a locally  stable  region. 

The  behavior  of  the  wave  disturbances  in  our  study 
resemble  certain  aspects  of  the  waves  observed  by 
Kri shnamurti  (1971a,  b)  along  the  South  Asia  easterly  jet 
during  the  Northern  Hemisphere  summer,  and  the  waves  that 
were  simulated  by  Colton  (1973)  in  a nonlinear  model.  The 
wave  disturbances  entering  upstream  of  the  easterly  jet 
regime  grow  and  increase  their  speed  and  wavelength  as  they 
approach  the  jet  maximum.  They  then  achieve  their  maximum 
amplitude  a considerable  distance  downstream  of  the  jet 
maximum  and  subsequently  decay  as  they  leave  the  strong  jet 


regime. 

The  simple  mechanistic  analytical  model  developed  in 
this  study  can  crudely  approximate  the  adjustment  of  tn 
moving  disturbance  to  the  local  stability  of  a var  r 
flow  within  the  dyriarn  it  a 1 ly  un  Suable  icy  On.  It 
that  this  simple  approach  could  be  used  t o 
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Although  this  study  employs  linear  equations,  it  is 
possible  to  estimate  qualitatively  the  nonlinear  effect  of 
the  waves  on  the  mean  flow.  Since  the  Reynold's  stress  term 
is  the  dominant  term  in  the  energetics,  most  of  the  mean 
kinetic  energy  is  removed  by  the  waves  on  the  downwind  side 
of  the  jet  maximum  between  where  the  wave  amplitude  is 
maximum  and  where  the  phase  tilt  of  the  waves  is  strongest. 
Clearly  these  synoptic  scale  waves  will  affect  the  amplitude 
and  phase  of  the  quasi-stat ionary  planetary  waves  which 
combine  to  form  the  easterly  jet  over  South  Asia.  These 
effects  can  be  studied  by  including  nonlinear  effects  in 
the  present  model.  It  may  also  be  interesting  to  extend  the 
present  model  to  a 3-dimensional  model  to  study  the  effect 
of  vertical  shear  and  ba roc  1 i n i ci ty  . 


A 
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APPENDIX 

Subroutine  POISDD 


Subroutine  POISDD  is  a direct  method  solver  of  the 
finite  difference  approximati on  to  Poisson's  equation 
(Sweet , 1971): 


v u(x,y)  = - f ( x,  y ) 


for  a<x<b 


c<y<d 


and  with  Dirichlet  boundary  conditions 
u(x,c)  = g1 (x) 

u(x,d)  = g£(x)  a<_x<_b  (2) 

u(a,y)  = g3(y)  c<y<b 

u(b,y)  = g4(y) 

We  define  a grid  on  the  rectangle  {(x,y):  a<_x<b, 
c<_y<_d}  by  selecting  two  positive  integers  M and  N,  such  that 
aX  = and  AY  = ^ (3) 

and  where  N must  be  a power  of  two,  and  the  points 


x . = a + ( i - 1 ) a x 

yj  = c + ( j-i ) ay 


1=1.  2, 


, M+l 


j = l,  2, N+ 1 , 


where  i and  j refer  to  the  x and  y grid  points,  respectively 
The  Poisson  equation  that  is  solved  by  POISDD  is  given 
by  Eq.  (2.10): 


?2  ft”  = ‘ ) - j(*'  ,V2$)  - Syij-—  - Df72V 
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with  boundary  conditions  Eqs.  (2.11),  (2.14)  and  (2.16).  The 
finite  difference  form  of  these  equations  are  given  by  Eqs. 
(4.10),  (4.14)  and  (4.17).  Here  T . t . is  the  finite  differ- 
ence  representation  for  the  tendency,  3’jj  /at,  and  is  the 
dependent  variable  in  our  Poisson  equation.  For  our 
Numerical  Model,  M = 1 0 6 and  N = 32;  t is  the  current  time, 
where  t = 4At  i = 1,2....  and  At  in  the  time  step. 

Our  system  of  equations  meets  the  requirements  of 
subroutine  POISDD,  except  for  the  outflow  boundary  condi- 

[ 

tion,  Eq.  (2.16). 


3 / 3 1 \ _ 3 / 3 ij;  1 \ 

3 1 ' 3 1 ' ~ Cr  3X^  3t  ' 

0 


(6) 


This  equation  is  the  Sommerfeld  radiation  condition  (Pearson, 


1974).  Clearly,  this  is  not  the  usual  form  of  Dirichlet 
boundary  condition.  Therefore,  in  order  to  use  POISDD  in 
our  numerical  model,  this  subroutine  had  to  be  modified 


incorporating  the  radiation  outflow  condition.  -as 

done  by  writing  Eq.  (6)  in  finite  difference  form: 


where  the  finite  difference  symbols  have  their  usual  meaning 
in  numerical  weather  prediction.  Subscript  i=l  refers  to 
the  outflow  boundary,  and  subscript  i=2  refers  to  the 
column  of  grid  points  next  to  this  boundary.  We  observe 
that  this  scheme  is  backward  in  time  with  upwind  space 


differencing.  Equation  (7)  is  conveniently  rewritten  such 
that 


Tt  - A t 

1 ,j 


T ^ 

2 ,.j 


T 1 , j ~ 1 + C " At/AX  + 1 + AX/ic  At)  ^ 8 ^ 

ro  0 

Since  all  the  interior  points  of  the  domain  must  satisfy  the 
non-homogeneous  Poisson  equation,  we,  then,  use  this  constraint 
to.  couple  the  radiation  boundary  condition  [Eq.  (8)]  with  the 


interior  points,  as  follows: 


T3J  - 2T2!j  * Tl!,i  + T2J  + 1 


tz!-+i  ' 2Tz!.i  * T2(.i-1 


'f2,r 


The  term  f.  . is  the  finite  difference  approximation  to  the 

I j J 

right  side  or  non-homogeneous  part  of  Eq.  (5)  or  Eq.  (2.10). 

2 

We  multiply  Eq . (9)  by  -(Ay)  , and  we  obtain 

T,^ . 21,*.  T , t . f t 

^ ^ + i + 2T7t.  - T?ti  1 

s s s 2.J  + 1 2 , j 2.J-1 


f2,j^ 


2 2 

where  s = ax  /Ay  . 


We  now  combine  Eq . (8)  and  Eq . (10),  and  we  get 


2 , j + 1 


+ 2 + r - 


V* 


t . OJ.  + 


T t-At 

1 Kj 


T2  , j + 1 _ f2,jAy  + s L 1 + c ^ J At/Ax. 

o 
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By  incorporating  Eq . (11)  in  subroutine  POISDD,  we 
have,  thus,  coupled  the  Sommerfeld  radiation  boundary  condi- 
tion to  the  interior  domain.  Also,  by  Eq.  (8),  we  now  have 
a complete  set  of  Dirichlet  boundary  conditions.  The  finite 
difference  Eqs.  (4.10),  (4.14)  and  (4.17),  thus,  form  a 
linear  system  of  equations  of  dimension  (M-l)  x (N-l ) for 
the  unknowns  T.*.:  2<_i<_M,  2<.j<N.  This  system  is  solved  in 

* J 

POISDD  by  the  Buneman  variant  of  the  CORF  (cyclic  odd-even 

reduction  and  factorization)  algorithm.  Buzbee  e_t  aj_.  ( 1 970  ) 

gtves  a complete  mathematical  description  of  the  algorithm. 

Sweet  (1971)  gives  a brief  description  of  the  subroutine. 

A function  of  two  variables,  B.  is  used  in  the 

' * J 

numerical  model  code  to  represent  a two-dimensional  array 

which  provides,  on  input  to  POISDD,  values  of  the  function 

f.  . as  well  as  the  specified  boundary  conditions.  On  output 
"l  * J 

from  POISDD,  it  provides  the  values  of  the  approximation 
T.C..  For  our  numerical  model,  we  approximate  the  outflow 

' * J 

boundary  condition  on  input  to  POISDD  as  follows: 


(12) 


On  output  from  POISDD,  we  correct  Eq . (12)  by 


(13) 
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